Finite-temperature phase transitions in quasi-two-dimensional spin-1 Bose gases 
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Recently, the Berezinskii-Kosterlitz-Thouless transition was found to be mediated by half-quantum vortices 
(HQVs) in two-dimensional (2D) antiferromagnetic Bose gases [Phys. Rev. Lett. 97, 120406 (2006)]. We study 
the thermal activation of HQVs in the experimentally relevant trapped quasi-2D system and find a crossover 
temperature at which free HQVs proliferate at the center of the trap. Above the crossover temperature, we 
observe transitions corresponding to the onset of a coherent condensate and a quasicondensate, and discuss the 
absence of a fragmented condensate. 

PACS numbers: 03.75.Lm,05.30.Jp,03.75.Mn,64.70.Tg 



I. INTRODUCTION 

The dimension of the underlying space has a profound im- 
pact on the existence of long-range order and phase transi- 
tions in a given system. In two-dimensional (2D) systems the 
long-range order and spontaneous symmetry breaking are for- 
bidden llll-yl]; however, 2D systems can exhibit a quasi-long- 
range order with algebraically decaying correlations [|4)-|7|]. 
The disordered high-temperature phase and the algebraically 
ordered low-temperature phase are separated by a topologi- 
cal phase transition corresponding to the unbinding of pairs of 
vortices and antivortices. This phase transition is referred to 
as the Berezinskii-Kosterlitz-Thouless (BKT) transition 
However, experimentally relevant examples often display ad- 
ditional features due to the finite-size effects [8], and in the 
trapped ultracold atomic gases where the BKT transition has 
recently been studied 191 4 1311 . the inhomogeneous density of 
the gas renders the superfluid state and the coherence proper- 
ties qualitatively different from those of the bulk systems Ill3[ - 

m. 

Spinor Bose gases lfl6l - [T8ll are especially interesting as they 
can in principle combine magnetic ordering, formation of a 
condensed component, and superfluidity. Due to the interplay 
of these competing orders, the antiferromagnetic spin-1 Bose 
gas is expected to host various exotic phenomena such as frag- 
mented condensates [19] and fractionalized topological ob- 
jects liolElll that are usually absent in the single-component 
systems. For example, a half-quantum vortex (HQ V) confined 
to a spin defect occurs in spin nematic condensates |20, 2^ 
and it has recently been created using Raman-detuned laser 
pulses ll23ll . In homogeneous 2D optical lattices, prolifera- 
tion of HQVs in spin-1 Bose systems due to thermal fluctua- 
tions has been predicted |24, 25], and the superfluid transition 
in two dimensions was found to be mediated by HQVs 1.26.] . 
In ferromagnetic spinor condensates, HQVs have shown to 
give rise to intriguing magnetization dynamics [tT]. Frac- 
tional vortices and the related BKT transitions have also been 
discussed in the context of '^He ll28ll and different nonconven- 
tional superconductors [|29l - [3Tll . Recently, HQVs have been 



observed in exciton-polariton condensates 03211 . 

While the connection between superfluidity and Bose- 
Einstein condensation is relatively well understood in the 
single-component Bose systems flHISHlSlj, the existence of 
spin degree of freedom in spinor Bose gases renders the rela- 
tion between superfluidity and long-range order more com- 
plicated and far less studied. In particular, the existence 
and the nature of the possible condensed component is not 
yet known in two dimensions. The recent experimental in- 
terest in spinor Bose gases with antiferromagnetic interac- 
tions (3?, 34*], advances in the evaporative cooling of opti- 
cally trapped atoms |35], and the nondestructive imaging of 
the local magnetization of spin-1 Bose gases ll36l IjtII sug- 
gest that the experimental realization of the finite-temperature 
phase transitions in quasi-2D spinor Bose gases may be pos- 
sible in the near future. Hence, we study the activation of dif- 
ferent topological defects associated with the superfluid tran- 
sition and determine the different degenerate components of 
quasi-2D antiferromagnetic spin-1 Bose gases. Our approach 
is valid in the regime where the thermal fluctuations are domi- 
nant and our results suggest that in this region, the condensate 
state is nonfragmented. 



II. FORMALISM 

To study the behavior of a spinor Bose gas near the crit- 
ical region we use a classical field (c field) to describe the 
highly occupied low-energy modes and a quantum field for 
the thermal modes with low occupation |38]. Previously, this 
approach has been successfully applied in studies of the BKT 
transition in scalar Bose gases [3 [lH Hit] as well as to 
predict other properties of dilute scalar Bose gases [40-4^. 
The coherence properties of spinor Bose condensates at finite 
temperatures have recently been studied using an alternative 
formulation of the c-field method 14711 . 



The dynamics of the c field can be described by the pro- 
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jected Gross-Pitaevskii equation (PGPE) 



-V{c^\^a\ 



(1) 

which generalizes the usual spinor Gross-Pitaevskii equa- 
tion (3 [I3l- In the PGPE, denotes a vector of spin-1 
matrices and V implements a projection into the subspace 
of the classical modes ifssll . The c field in Eq. ([T]) is writ- 
ten in the basis consisting of the Zeeman substates such that 
^ c — (V'c a)' Q; = 1, 0, —1. The single-particle operator 
is given by 
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The coupling constants Cq and C2 are given by Cq = Anh^ (ao+ 
2a2)/3mandc2 = 47r?i^(a2~aQ)/3m,wheremis the atomic 
mass and Cq and 02 are the s-wave scattering lengths in the 
total hyperfine spin channels F = Q and F = 2, respec- 
tively 1 17]. Antiferromagnetic interactions imply C2 > 0, and 



we take Oq — 46a^, 



I2 = 520^, and m = 3.86 x lO^^e j^g^ 



according to ^■^Na 1 17]. The Bohr radius is denoted by a^. 

In the quasi-2D limit, lo^ <C and we choose = 200 x 
Wj^. Harmonic oscillator lengths in axial and transverse direc- 
tions are denoted by = ^ h/mujz and Cj^ = ^Jh/mujZ. 
The scatterin g ca n be treated as three-dimensional as long as 
Oq, 02 <C ll48l l49tl . This condition is satisfied with the pre- 
vious choice of and uj^ if we take = 27r x 10 Hz which 
is in the realm of the current experiments. Since we consider 
a quasi-2D situation, the PGPE and all c fields are defined in 
a 3D space. The spatial vector of the 3D space is denoted by 
r and is a 2D vector in the x — y plane. Summation over 
repeated indices is implied. 

In the PGPE, the c-field region C is defined by the energy 
cutoff e^yj such that C = {n | e„ < Ecutl' corresponding to 
the spectrum of the single-particle operator /iq. The c fields in 
Eq.O can be expressed in terms of the eigenstates of /iq 



(2) 
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The PGPE corresponds to a microcanonical system in which 
the stationary probability distributions are determined by the 
total energy of system, and the temperature and the chemi- 
cal potential are computed as ensemble averages. We use the 
ergodic hypothesis to replace all ensemble averages with the 
corresponding time averages. Using the ergodic hypothesis, 
thermodynamical quantities such as the temperature and the 
chemical potential can be computed dynamically [38, 41, 50]. 

Let us briefly discuss how to generalize the single- 
component calculation of the temperature and the chemical 
potential |38, 41, 50] to the spin-1 case. The PGPE dTJ arises 
from the Hamiltonian 

Ha=J dr [^iK^^ + + ^-li^lT^^f] (3) 

for which the canonical coordinates can be defined such that 
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Qa,n ^ AT ('^a.n^"'^Q,7^) P^^ 



(4) 



where „ are the coefficients in Eq. (|2|i. The canonical co- 
ordinates are collectively denoted by F = {Q^ „}. Ac- 
cording to a general theorem ifsill . the temperature can be cal- 
culated as 



(5) 



where the first identity is the standard definition of the tem- 
perature of a microcanonical system. The derivative operator 
D = {e„9/9p }, determined by the coefficients {e„}, and 
the vector field Xj, can be chosen freely as long as they sat- 
isfy the conditions SIHIll 



VHc ■ Xj, = l and VN^ 



(6) 



where Nc = /dr I'I'cP is the total number of the c-field 
atoms. The vector field Xj, satisfying the above constraints is 
given by [41., .50.] 



(7) 



withAjv = {VNc ■ VHc)/\VNc\'^. Straightforward choices 
for the vector operator T) are Dp = {0, 9p } and = 
{5q , 0} jfsO]. The temperature is independent of the choice 
the derivative V, and the two different choices serve also as 
a check for the numerical implementation. The average in 
Eq. (|5]l is computed as a corresponding time average. 

The preceding formulation can be applied when the only 
conserved quantity is the total particle number Nc- In the 
present case, also the angular momentum is conserved. Trans- 
forming to the coordinate system with zero total angular mo- 
mentum, the angular momentum conservation does not ap- 
pear in Eqs. ^ and Q ifsol I52I] . For spinor Bose gases, 
the conservation of the total magnetization also needs to be 
taken into account. In the antiferromagnetic case, however, 
the total magnetization is zero and can be neglected in light 
of the previous argument. Using the definition /i/fcgT = 
~~{dS / dN) ^, also the chemical potential can be computed 
by interchanging the roles of He and Nc- Computationally 
efficient formulation for the different terms in Eqs. (|5]l - (|7]i 
proceeds in an analogous way to Ref . llsoll . 

The number of atoms outside the c-field region can be com- 
puted self-consistently using the Hartree-Fock-Popov (HEP) 
approximation il [H m m. The full field operator con- 
taining the c-field part -0^ ^ and the incoherent part 5(t>i,a is 



denotedby $^ = V'c,q+<5' 



'La 



We assume that terms such as 



and all their complex conjugates 
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vanish. This leads to the HFP single-particle energies 11531 



e^(fc, r) 

eo(fc,T-) = 
e_(fc,r) 

where n„ 
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+ VtAr) - H + Coin + no) + c^in^ + n_), 

(8b) 

+ Kr(^) - /i + Co(n + n_) 

+ C2(2n_ +no-n_^), 
(8c) 

= l^cal' + and n = (fi) = 



+ + n_. The occupation number n^^' = 
for the incoherent atoms can be computed from the Bose- 
Einstein distribution using the semiclassical integral in a 3D 
phase space 



dfc 



(27r)3 f,eJk,r)/k^T _ i 



(9) 



The quasi-2D nature of the system is taken into account by 
treating the axial modes discretely in the semiclassical inte- 
gral U41 . The energy cutoff e^^^ introduces a spatially depen- 
dent low-energy cutoff to the semiclassical integral (|9]l (see 
Ref. Q). 



III. TOPOLOGICAL DEFECTS AND NEMATIC ORDER 

The c field in Eq. ([T]i is written in the z-quantized basis 
^ = (V'q), a = ~1, 0, 1, but the nematic properties of anti- 
ferromagnetic Bose gases are more conveniently expressed in 
the Cartesian representation lfT6ll54ll ^ = (^Aa), a = x,y,z. 
The transformation is given by tp^ = {ifi-^ — tp_^)/^/2, 
i'y = + 4'-i)/V^, and -0^ = and the nematic or- 
der is described by the spin quadrupole moment 1,54] = 

(■0*^6 + i'b'^a)l{'^\^?)- In general, Q'") has three distinct 
nonzero eigenvalues and the local magnetic axis fi is defined 
as the eigenvector associated with the largest eigenvalue. For 
T/;^ = 0, the magnetic axis is confined into the x — y plane and 
we refer to such case as the "in-plane" nematic. 

In the polar phase which corresponds to identically vanish- 
ing local spin iflTll . the Cartesian representation gives ^ ~ 
Y^e'^n and the HQV corresponds to a defect where both 9 
and n have a tt winding about the core of the defect. Further- 
more, the polar phase allows also the existence of skyrmions 
which have finite energy and are characterized by the second 
homotopy group of the order parameter space. In the Ap- 
pendix, explicit expressions for the HQVs and skyrmions are 
presented. The polar phase has a local Z9 invariance corre- 
sponding to {e,h) (9 + TT, -n) imiilEll. This im- 
plies that defects with opposite topological charges cannot be 
distinguished and therefore we define the sign of the HQV 
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Figure 1: (Color online) Instantaneous density of the c-field atoms 
corresponding to the out-of-plane nematic phase. Half-quantum vor- 
tices and antivortices are denoted by black + and — symbols, re- 
spectively. Skyrmions are marked with black dots. The black line 
denotes the boundary outside which Uc/fitot < 0.1, where is the 
condensate density, and ntot is the average total density of the c-field 
atoms. The thermal wavelength is denoted by A. In the instantaneous 
density, the 2-dependence is integrated out. The temperature is given 
with respect to the critical temperature of the corresponding quasi- 
2D ideal Bose gas (see Section III). 



from the polarization of the vortex core 112211 . Furthermore, we 
do not di sting uish between skyrmions with opposite winding 
numbers Ipql . An example of the thermally activated HQVs 
and skyrmions is shown in Fig. [T] where the instantaneous z- 
integrated densities of the c-field atoms at different tempera- 
tures are depicted. 

We consider two phases, the in-plane nematic with n^ = 
n^ _i and n^, g = and an "out-of-plane" nematic with 
^c+i ~ ^c,o ~ "^c -1 = 0-33 ± 0.06. Here refers to 
the average number of c-field atoms in the component a, di- 
vided by the total number of atoms Nc in the c-field region. 
The average numbers ^ corresponding to the different data 
points in Figs. |2] and [3] vary between the aforementioned lim- 
its. In both cases, we take Nc — 15000 and choose the energy 
cutoff as e^jjj = 126 hw^ i^cnt = 122 hw^) for the in-plane 
(out-of-plane) nematic. For these values of the cutoff energy, 
only the lowest axial mode becomes populated and the c-field 
region is in the ground state with respect to motion in the z 
direction. Nevertheless, the incoherent region atoms typically 
occupy several axial modes. 

As indicated in Ref. |25], the in-plane nematic phase arises 
in the spin-1 case as a result of a large negative quadratic Zee- 
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man shift (for a discussion how the negative shift is physically 
realized, see Refs. iIzsIIstIi ). In the case of the PGPE, elimina- 
tion of the a = component corresponds to leaving the a — 
component empty in the initial state. The quadratic Zeeman 
shift can be absorbed in the single-particle energies since it is 
the same constant for the a = ±1 components. This allows us 
to the treat the in-plane and out-of-plane cases at equal foot- 
ing, assuming only that the Zeeman shift is large enough to 
eliminate the a = component at the relevant temperatures. 

The ensemble averages are calculated as corresponding 
time averages such that the system is allowed to thermalize 
for period 50 x 2tt/ujj_ and the time average is computed 
from 1250 equally spaced samples. The sampling interval is 
50 X 2-K juj ^ (100 X 2tt/uj^) for the in-plane (out-of-plane) 
nematic phase. The randomized initial states are taken from 
the polar phase corresponding to zero magnetization. Oth- 
erwise the numerical implementation follows the description 
of Refs. issils^ . In the HFP calculation for the in-plane ne- 
matic, we assume that there are no thermal atoms in the a ~ 
component. 

We keep the cutoff energy fixed, which causes the total 
number of atoms iV^^j to increase with increasing tempera- 
ture. To accommodate to the varying atom number, we scale 
the temperature by the critical temperature TJ, of a quasi-2D 
ideal Bose gas corresponding to the same total particle num- 
ber 1 15]. For the in-plane (out-of-plane) nematic phase, Tg 
corresponds to the critical temperature of two (three) inde- 
pendent ideal Bose gases. For the ideal Bose gas, there is a 
standard relation between the temperature and the occupation 
number of the thermal component 1 15. .5911 

where are the 3D harmonic oscillator energies. We calcu- 
late Tq numerically by solving the equation N{Tq) — 7iVtot 
where 7 = 1/2 for the in-plane phase and 7 = 1/3 for the 
out-of-plane phase. We denote the reduced temperatures by 
T' and the bare temperatures by T. In experiments, a pre- 
cise control of the total atom number at different temperatures 
is difficult and hence the approach presented here is likely to 
describe well the realistic experimental conditions. 



IV. CONDENSATE AND QUASICONDENSATE 

The existence and the nature of the condensate and the qua- 
sicondensate components in antiferromagnetic Bose gases are 
particularly interesting due to the possibility of a fragmented 
condensate at zero temperature [191. Since the fragmenta- 
tion in this case corresponds to the condensation of composite 
bosons to the |fe = 0) state in the momentum space, it seems 
that also the fragmented condensate is destroyed by the ther- 
mal fluctuations in a homogeneous 2D system. In addition, 
the thermally activated HQVs render the single-mode approx- 
imation used in Refs. lfl9[ l60ll invalid and it is a nontrivial 
question whether the fragmented condensate can exist in 2D 
at finite temperatures. In this work, the presence of a signif- 
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Figure 2: (Color online) Condensate fraction Ng/Ntot for the in- 
plane nematic (blue circles) and the out-of-plane nematic (red trian- 
gles) phases as a function of the reduced temperature. The quasicon- 
densate fraction is almost identical for the in-plane nematic (green 
squares) and out-of-plane nematic (black asterisks) phases. Note that 
the definition of the reduced temperature differs between the in-plane 
and the out-of-plane phases (see text for details). (Inset) The conden- 
sate fraction near the temperature corresponding to the onset of the 
condensate. The dashed line is a guide for the eye. 

icant thermal component renders a direct comparison to the 
zero temperature single mode calculations difficult. 

In the homogeneous 2D case, algebraic order is expected in 
the paired state coiTesponding to 9 = $o^o ~ '^^+1^-1 
and inspection of the correlation function (G^(r^)0(r^)) 
could shed light on the superfluid properties of the antifer- 
romagnetic spin-1 Bose gases. In this work, we are interested 
in the existence and the nature of a condensed component in 
spin-1 superfluids and consider therefore the one-body density 
matrix p^^\ra;r' (3) = (^^(f which can be sam- 
pled using the time averaging. Under the previous assump- 
tions it separates into two parts containing the c-field part and 
the incoherent part. At low temperatures, we find that p'^^^ 
has only a single macroscopic eigenvalue A^q and we refer to 
-^o/^tot by the generic name "condensate fraction." 

Above the critical temperature of condensation, p^^^ has 
several large eigenvalues although their fraction of N^_^^ be- 
comes vanishingly small. This thermally induced fragmenta- 
tion ll60ll is, however, different from the fragmentation due to 
the ordering in the spin sector. Our results seem to be con- 
sistent with the idea of a hierarchy of transition temperatures 
such that the formation of a coherent condensate is followed 
by ordering in the spin sector, leading potentially to a frag- 
mented condensate in the zero-temperature limit |60]. The 
condensate fraction is shown in Fig.|2]as a function of the re- 
duced temperature T' = T/Tq. 

For the scalar Bose gas, the quasicondensate component 
can be defined using the correlation function C{rj_) = 
2(0t(r^)0(r^))2 _ ([^t(r^)0(r^)]2) |i5[|6l, describing 
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ponent tends to be localized to the region where HQVs and 
skyrmions are rare (see Fig.[TJ. 

Since the reduced temperature for the onset of a condensate 
is the same for both nematic phases within the numerical accu- 
racy, it is natural to ask if it is caused by the condensate deple- 
tion due to the incoherent atoms. Although the number of the 
incoherent region atoms is large near the onset of the conden- 
sate, the same onset temperature for the condensate is found 
if only the c-field atoms are taken into account. In Fig. [3] the 
fraction Nq/Nc is shown as a function of the reduced tem- 
perature, indicating that the onset of a large eigenvalue takes 
place at equal temperatures for both nematic phases. Hence, 
the onset of a nonzero condensate fraction at the same tem- 
perature for both nematic phases does not depend on the de- 
pletion of the condensate due to the incoherent region atoms. 



V. PROLIFERATION OF TOPOLOGICAL DEFECTS 



Figure 3: (Color online) The largest eigenvalue of the one-body den- 
sity matrix (A^o) normalized by the number of the c-field atoms (Nc ) 
as a function of the reduced temperature. The inset shows the same 
quantity zoomed to the temperatures corresponding to the onset of a 
large eigenvalue. The out-of-plane nematic is denoted by (red) tri- 
angles, and (blue) circles correspond to the in-plane nematic phase. 
The dashed line is a guide for the eye. 

the part of the system with reduced total density fluctua- 
tions. In the spinor case, we may consider an equivalent quan- 

tity CK) = 2{^Ur^)^air^r ~ (['^J.K)$„K)]^> 
where the z dependence in the expectation values is inte- 
grated out. In the spin-1 case analyzed here, we observe 
within the HFP approximation that the fraction of the com- 
ponent with suppressed total density fluctuations given by 
/ dr^ ■\/C(rj_)/A^jQ4, decreases very slowly with increasing 
temperature. This illustrates the role of the thermally induced 
intercomponent density fluctuations. 

We determine the quasicondensate component by consider- 
ing the total density fluctuations restricted to the c-field region 
and define the quasicondensate density as 

n,,(r) = V'2(l*o(r-)P>2-(l*oWI^>. dD 

For the parameters discussed in Section |III1 the c field is of 

the form ^cl*") ~ '^c{'''±)fo{^)' where (Pq{z) is the lowest 
harmonic oscillator state. Hence, the z dependence in n^^ir) 
can be integrated out. The quasicondensate fraction is shown 
in Fig. |2]and it persists at the temperatures where the conden- 
sate fraction becomes negligible. The critical temperature for 
the formation of the coherent condensate as well as the onset 
of the quasicondensate are the same for the in-plane and the 
out-of-plane nematic phases, and they take place at tempera- 
tures T = 0.97 ± 0.02 and T' = 1.16 ± 0.04, respectively 
(from Fig.|2]i. In addition, the quasicondensate fraction is es- 
sentially the same at equal temperatures in both cases. Due to 
the reduced total density fluctuations at all temperatures, the 
quasi-condensate component is delocalized to the entire spa- 
tial extent of the c-field atoms whereas the condensate com- 



In trapped atomic gases, the characteristic feature of the 
crossover from a BKT type of superfluid to a normal fluid is 
a proliferation and a subsequent propagation of free vortices 
from the edge of the cloud to the central region of the trap. 
Since HQVs are nonsingular defects, they persist at the edge 
of the cloud to relatively low temperatures (Fig. [Hi and the 
system can be considered to have concentric shefls of normal 
fluid and BKT superfluid with the center of the trap occupied 
by the condensate. We analyze the BKT crossover by studying 
the HQV occupation probability density 13. An estimate 
for the crossover temperature is obtained from the temperature 
at which becomes nonzero near the center of the trap (see 
Fig. |4]i. From this analysis, the BKT crossover takes place 
roughly at the reduced temperature T^^^^ = 0.82 ± 0.05 for 
the in-plane nematic and at Tg,^^ = 0.89 ± 0.04 for the out- 
of-plane nematic phase. 

Since the reduced BKT crossover temperature obtained pre- 
viously has a rather large uncertainty, we cannot conclusively 
determine the relation between the crossover temperatures for 
the two nematic phases. However, since the order parame- 
ter has a different symmetry in the in-plane and out-of-plane 
phases, there is no reason to assume that the crossover tem- 
perature is the same. For the in-plane nematic phase the sym- 
metry is reduced to [?7(1) x S^]/TL2 while in the case of an 
out-of-plane nematic it is \U{1) x S'^]/K2, allowing the ex- 
istence of skyrmions which in the homogeneous case render 
the system spin-disordered. In a finite-size system, the ther- 
mal activation of skyrmions depends on the characteristic size 
of skyrmions compared to that of the system, and we find that 
skyrmions start to appear only at relatively high temperatures 
near the BKT crossover, see Fig.|4] 

The effect of skyrmions to the crossover temperature Tbkt 
can be illustrated by considering the statistical probability for 
the activation of a skyrmion or a pair of HQVs. The probabil- 
ity is proportional to the Boltzmann factor cxp(— AJ^/fcgT) 
where AF — A£ — TAS" is the free-energy change associ- 
ated with the creation of a given defect. The critical temper- 
ature for the activation of different defects can be estimated 
from the condition AJ^ = 0. In a uniform 2D system, the en- 
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Figure 4: (Color online) Radial probability density for detecting (a) a HQV in the in-plane nematic, (b) a HQV in the out-of-plane ne- 
matic, and (c) a skyrmion in the out-of-plane nematic phase. The reduced temperatures are given by (from top to bottom in each panel) 
1.05, 0.95, 0.82, 0.63 in (a), 1.00, 0.94, 0.89, 0.64 in (b), and 1.00, 0.96, 0.92, 0.84 in (c). 



tropy change associated with the creation of a skyrmion can be 
approximately evaluated as 10,0- IHl ^S^^ — 2fcgln(€/ro) 
where £ is the linear size of the system and Tq the characteris- 
tic size of a skyrmion. The skyrmion energy E^^. — Anh^g/m 
is finite and independent of the size of the skyrmion (see Ap- 
pendix). Hence, the free energy is always negative for a large 
enough system and skyrmions exist at all temperatures in the 
thermodynamical limit |63]. 

Since the BKT transition is driven by the dissociation of 
bound pairs of vortices and antivortices ll64ll . we estimate the 
critical temperature using the Boltzmann probability for the 
appearance of free vortices [6, 7] . In the Appendix, the en- 
ergy of a free HQV is shown to be ^nqv = '^'^i^l^i)^ 
where is the size of the HQV core. The corresponding en- 
tropy change is AS'jjgy = 2fcg ln(7]£/rj^), where rj < 1 m 
the presence of skyrmions due to screening. This would re- 
sult in a higher critical temperature for the activation of free 
HQVs in the out-of-plane nematic phase where skyrmions are 
allowed. In nonuniform finite-size systems, the mechanism 
is different since the skyrmions appear only near Tbkt- The 
thermal fluctuations generate skyrmions first at the boundary 
of the cloud and since the skyrmion energy is independent of 
its size, this process is not strongly affected by the pre-existing 
HQVs. The generation of skyrmions can prevent the thermal 
fluctuations from breaking the HQV-anti-HQV pairs, thereby 
giving rise to the higher crossover temperature. 

The crossover temperature can also be studied using the 
2D phase-space density nJ,^°'A^, where n*,^"' is the aver- 
age 2D total density of the c-field atoms at = and 
A = ^J^jrWJmk^ . We observe that the phase-space density 
takes roughly the value 25 for both nematic phases at the re- 
spective reduced crossover temperatures (see Fig.|5]l. This re- 
sult is to be contrasted with the single-component case where 
the transition to the superfluid phase takes place when the 
phase-space density is larger than the critical value 

n<:rt'A'-ln(C/g), (12) 

where g = ^/Stt a/a^, a is the s-wave scattering length in 
the single component system, and C ~ 380 ifTst l6lll . For 



the spin-1 system considered here, the corresponding 2D cou- 
pling constants can be defined as Cq — \/87r(ag + 2a2)/3az 
and C2 = V8TT{a2 — aQ)/3az- The parameters defined in 
Sec. in yield 5o = 0.028 and ^2 = 0.0011. The value 
of Cq is comparable to the experimental value g = 0.02 in 
the single-component system of Ref. |13] where the spin- 
dependent interaction is absent. A simple-minded application 
of the scalar case condition (fT2l) using either of the coupling 
constants Cg and with C ~ 380 yields much lower values 
than n^^°*A^ = 25. This suggests that if a condition analo- 
gous to Eq. (fTZt exists for the spinor case, its form is different 
from (inii. 
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Figure 5: (Color online) The phase-space density n!,^°'A^ as func- 
tion of the reduced temperature for the in-plane nematic (blue circles) 
and out-of-plane nematic (red triangles). The solid lines denote the 
BKT crossover temperatures Tbkt = 0.82 and Tbkt = 0.89 for 
the in-plane and the out-of-plane nematic phases, respectively. The 
dashed line corresponds to the phase-space density n|.J^'A^ = 25. 

An important check for the observed crossover tempera- 
tures is the superfluid density which is predicted to change 
in the spinor case noncontinuously across the phase transition 
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by the amount |12 



(13) 



Hence, the universal jump in the superfluid density is four 
times larger than in the single-component case. To use this 
property to check the consistency of the crossover tempera- 
tures, an independent computation of the superfluid density is 
required. Since the system is inhomogeneous, the central part 
of the system is typically in the superfluid state while the outer 
part corresponds to normal fluid. This renders methods such 
as the computation of the helicity modulus ll26l l65ll inappli- 
cable since they require a uniform system without coexisting 
phases. 

The HQVs are nonsingular vortices and the nonclassical 
moment of inertia [12.1 does not capture the BKT crossover 
The phenomenological models for the trapped systems in the 
single-component case fl2l [T5ll make use of the condition (fTZt 
and assume explicitly a sudden change in the superfluid den- 
sity by the universal value 2m?k^T^^^ / {-kH^) . Hence, there 
is a clear incentive for further investigations of the superfluid 
properties of spinor Bose gases, in particular for the deter- 
mination of the superfluid fraction without making use of 
Eq. ( fT3] ). We note that the condition n^^°^\^ — 25 is con- 
sistent with Eq. ( fT3] l since the crossover temperature yields 
(Pg = nirig) n^/n'-^^^ = 0.64. Using the scalar condition ( fTST i 
with g equal to either Cq or C2 gives J^s/^crit' ^ 1' indicating 
that the scalar condition (fT2] l is not vaUd for the spinor case. 

We note that the uncertainty in the determination of the 
BKT crossover temperature in terms of the reduced temper- 
ature Tg^^/To allows in principle the crossover temperature 
to be the same for both nematic phases. However, if there 
is a critical value for the phase-space density independent of 
the type of the nematic ordering, then the data in Fig. |5]yield 
different reduced crossover temperatures if they are below the 
condensation temperature T' w 0.97. We also note that it 
is numerically difficult to distinguish between skyrmions and 
merons when there are large fluctuations in the direction of 
the magnetic axis n(r^), but in an analogy to the homoge- 
neous 2D situation, we refer to these out-of-plane defects as 
skyrmions. In the in-plane case, skyrmions are not allowed 
but, instead, integer vortices corresponding to winding 27r of 
in the magnetic axis around the vortex core can take place. 
Such vortices seem to remain suppressed suggesting that they 
are irrelevant for the BKT crossover 



VI. DISCUSSION 

We have analyzed the realization of the BKT transition 
in antiferromagnetic spin-1 Bose gases under typical exper- 
imental conditions. We have found a hierarchy of crossover 
temperatures corresponding to the onset of a quasicondensate 
at a high temperature and the formation of a coherent con- 
densate at a lower temperature, followed by a BKT-type of 
crossover to a superfluid state as the temperature decreases. 
The investigation of the probability density for the excitation 
of skyrmions and half-quantum vortices did not unequivocally 



determine the relation between nematic ordering and the mag- 
nitude of the crossover temperature. The subsequent inspec- 
tion of the 2D phase-space density at the center of the trap 
suggested that the crossover temperature expressed in terms 
of the reduced temperature could be slightly higher for the 
out-of-plane phase. Further investigations are still needed to 
confirm this scenario. The finite size of the system is mani- 
fested as a finite activation temperature for the skyrmions and 
the thermal fluctuations start to generate skyrmions only near 
the crossover temperature. It remains an open question if an- 
other crossover to a fragmented condensate takes place in the 
zero-temperature limit. 

We expect that the fractional population of different Zee- 
man sublevels can be controlled using rf pulses and magnetic- 
field gradients |66] to allow the experimental realization of 
the in-plane and the out-of-plane nematic phases. Using the 
time-of-flight imaging combined to the Stern-Gerlach separa- 
tion of the different Zeeman sublevels, formation of the con- 
densate component can be observed. The ferromagnetic cores 
of the HQVs could be detected by imaging the magnetization 
of the gas (37^, and the same technique can in principle 
be extended to image directly also the spin quadrupole or- 
der |36t ll^l- Interference experiments |10] can also be 
useful to demonstrate the existence of free vortices around the 
BKT crossover temperature. 
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Appendix: Energies of skyrmions and HQV-anti-HQV pairs in 
uniform systems 

The skyrmion configuration can be represented on a Carte- 
sian basis such that $ — y/pe'^n, with 



n = (sin cos sin (3{p) sin tp, cos (3{p)) , 



(A.l) 



where {p, ip) denote the polar coordinates and function /3{p) 
satisfies the boundary conditions (3{0) = and f3{p) = n 



for p > 



A meron (half-skyrmion) is obtained with 



/3{p) = -K 12 for p > Tq. We assume a uniform system such 
that the density is a constant for skyrmions and HQVs. For 
the skyrmion configuration (lA.ll i. the C/(l) phase Q can be 
taken to be a constant. 

The low-energy theory for the polar phase is the nonlinear 
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cr model (NLcrM) of the form ||2^ 

C^^ I dr^[{\7h)' + {V0)'], (A.2) 

where the superfluid stiffness is K = h^g/m. The NLcrM has 
a conformal invariance such that the energy of a skyrmion 
can become independent of the size Tq. Furthermore, all 
configurations of the form ( I A. II ) satisfy the condition £^y_ > 
AttK [63]. Hence we can take the energy of the skyrmion to 
be ^ AnK. 

Outside the vortex core for p > r^, the HQV configuration 
corresponds to 9 — if / 2 and 

n = (cos|,sin|,0). (A.3) 



Assuming that the systems has linear size £, substitution 
of (IA.3I I into (IA.2b yields the energy ^nqv = ^ ^^i^/^i)- 
The HQV energy does not include the contribution from the 
vortex core which is negligible in the thermodynamical limit. 
Furthermore, the usual arguments [59.1 can be used to con- 
clude that the energy of a HQV-anti-HQV pair is in the lead- 
ing order f pHQv — ""-^ ^'^('^/'"i)' where d is the distance be- 
tween the vortex cores. 
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